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ABSTRACT 



In the lattice formulation of the Heavy Quark Effective Theory (LHQET), the 
"classical velocity" v becomes renormalized. The origin of this renormalization 
is the reduction of Lorentz (or 0(4)) invariance to (hyper)cubic invariance. The 
renormalization is finite, depends on the form of the decretization of the reduced 
heavy quark Dirac equation, and can persist down to zero lattice spacing. For the 
Forward Time - Centered Space discretization, the renormalization is computed 
both perturbatively, to one loop, and non-perturbatively, using an ensemble of 
lattices with (3 = 6.1 provided by the Fermilab ACP-MAPS Collaboration. The 
estimates of the leading multiplicative shift agree reasonably well, and indicate 
that for small classical velocities, v is reduced by about 20-25%. 



I INTRODUCTION 

In the Isgur-Wise heavy quark limit 1 ' 2 , there are many relations between the decay 
constants and form factors of particles containing a heavy quark. Among the most striking of 
these conclusions from the heavy quark spin-flavor symmetry is the fact that in the Mq — ><>o 

limit, a single form factor, the Isgur-Wise universal function describes all semileptonic decays 

of one meson containing a heavy quark into another, such as the process B — >D * lv . From the 

first, it has been emphasized that although the heavy quark spin-flavor symmetry suffices to infer 
the existence of this function, its calculation requires non-perturbative techniques, such as lattice 
gauge theory 3 . Several such calculations have been carried out. A calculation by the present 
authors used a lattice implementation of the heavy quark effective theory 4 , and two other lattice 
calculations treat the heavy quarks as Wilson fermions with a small hopping constant, but do not 
directly implement the heavy quark limit ' . 

On the lattice 4 , as in the continuum 1 ' 2 ' 7 , the Isgur-Wise limit entails the introduction of 
a "classical velocity", v, normalized to 1, which appears in the decomposition of the momentum 
of a heavy particle and the reduced Dirac equation of the heavy quark field: 

P = Mv + p 

(1) 

-ivDh iv \x) = 

In the continuum, the velocity that appears in these two contexts is the same. However, on the 
lattice this is not the case. In this article we explain the origin of this new renormalization and 
describe two calculations of its magnitude. One calculation is a one-loop perturbative calculation. 
This follows the analysis of Aglietti on a variant lattice Dirac operator. The second calculation 
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is non-perturbative. It is based on a computation of the shift in the energy of a meson containing 
a heavy quark, as measured by the fall-off of its Euclidean space propagator, for a given shift 
in its residual momentum. 

The incorporation of the above M — ><>o limit into a Lagrangian is accomplished by 
factoring out a phase which is singular in the M — >°° limit and defining a reduced field 7 

1 + Y ' V h^(x) =e iMv - x 1 + Y ' V V (jc) (2) 
2 2 

In the limit, the Lagrangian for h^ v \x) becomes 

S£ (v) = h (v \x) ivD h {v \x) = h iv \x) iv[d - igA(x)] h (v \x) (3) 

The propagator of the reduced field h (v \x) is 

S (v \p) = L— — (4) 

vp - £ (v) (p) 

where in perturbation theory E ^ v \p) is the sum of all 1-particle irreducible self mass diagrams, 
as illustrated by Figure 1. 
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Figure 1 — 



The proper self mass shift of the reduced heavy quark propagator 



II RENORMALIZATION OF THE CLASSICAL VELOCITY 

If £ ^ v \p) is Taylor expanded in powers of the residual momentum 



(v),_x - 3E (V) 



L KV, (p) = m + 



(5) 

p=0 



the first term gives a mass shift. This is without physical significance and can be eliminated by 
a redefinition of the reduced heavy quark field by a non- singular phase which is independent of 



the heavy mass 9 . 



h {v \x) ^ e imv ' x h {v \x) ^ 



This corresponds to a finitely different breakup of the total 4-momentum 
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P = Mv + p — > Mv' + (p + mv) ( 7 ) 

In the heavy quark limit, v and v' become equal. The renormalized classical velocity is 
inferred from the small residual momentum behavior of the full heavy quark propagator 



S< v >(p) = ±— (8) 

- m + v (ren) -p + - 



This identifies the renormalized classical velocity as: 



v <™, . zM(V(i _ ^ 



X. = 



p=0 



The wave function renormalization constant Z*^ is determined by the requirement that the 
classical velocity has unit (Minkowski space) normalization both before and after the shift 



v {ren)1 = v 2 = Z^\v - X) 2 = - 1 d°) 



In the continuum, X must be proportional to the only available 4- vector, v . Therefore, 
it modifies the inverse of the heavy quark propagator to first order in p only by an overall 
multiplicative factor, which amounts to a wave function renormalization of h^ v \x), but gives 
no velocity shift. 



On the lattice, because of the reduced Lorentz or rotational symmetry, X need not be 
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proportional to . For example, it can have a term proportional to , which has the same 0(4) 



transformation properties as v , but is linearly independent of it. When the coefficient X is not 
simply proportional to v , the classical velocity is shifted. 

An equivalent procedure for calculating the renormalized classical velocity, but which is 
better suited to non-perturbative simulation, is to consider the rate of fall-off of the reduced 
propagator in Euclidean time for fixed spacial momentum. 

E^\p)t (11) 



S(t,p) ~ zMe-^'W 



Allowing for a possible residual mass, the propagator decay rate is 



E (v \p) = limV(M + m) 2 + ((M + m)v {phys) + p) 2 - M\j\ + v^ i} 

■ (12) 

(phys) 



V 



To help keep track of the context, we will use the phrase the physical classical velocity, v 



(phys) 



to refer to the renormalized classical velocity, when it is defined non-perturbatively through the 
rate of decay of a propagator in Euclidean time. 
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The logarithmic derivative of the asymptotic propagator 



dS^(t,p)/d Pi \ p = Q _ dE ^(p) 



S (v \t,p=0) 



v i 



(phys) 



p=0 



(phys) 



t = 



~(phys) 
i 



(13) 



directly extracts the physical classical velocity. 



(phys) 



(l,v 



(phys) \ 



1 



r (phys) 



(14) 



We now proceed to two calculations of the renormalized classical velocity The first uses 
one-loop perturbation theory, and the second uses a lattice simulation. In both, we use the 
forward time - centered space discretization of the heavy quark Dirac equation 4 



v [ U 4 (x,x + i)S (v \x + t,y)-S (v \x,y) ] 



+ £ [f/^,^ + A)5 (v) (^ + A,v) - [yx,*-A)5 (v) (*-A,y)] - 5(x,y) 



(15) 



III ONE-LOOP RENORMALIZATION OF THE CLASSICAL VELOCITY 

The one loop the proper self mass is given by the diagrams shown in Figure 2. The point 
interaction only gives rise to a residual mass, and so can be ignored. With an infrared regulator 
X, the second diagram gives 



Figure 2 — 



The 1-loop contribution to the heavy quark proper self mass 



£ (2) -8 2 C 2 



d 4 k 



(2n) 



-V^p- 2k)A (X) (k)S^(p - k)V^p-2k) 



(16) 



The integration domain is periodic, [- n /a ,n /a] , and the bare propagators and vertices are 



S ' (p) = \e 4 - 1/ + 2_ sin p. a 



a 



A a>w = -2 £ sm -y- + 



iv e 



l « a/2 , V- COS 



(17) 



The factor of the quadratic Casimir invariant, C 2 , which is (N - \)I2N for SU(A0 color with 
quarks in the fundamental representation, arises from the color contraction in the loop integral. 

The Taylor series coefficient of p 1 is a non-linear function of v and is furthermore 



asymmetric between space and time (because of the asymmetry of the lattice action, which is a 



centered difference in space but a forward difference in time). Thus there is a classical velocity 

shift as well as a multiplicative renormalization constant Z^. 

As has been noted, the symmetrical first difference results in fermion doubling. However, 
unlike in the usual fermion situation, the contributions of the secondary modes automatically 
vanish in the zero spacing limit. The reason for this is that each is accompanied by an undoubled 
gluon at the edge of the Brillouin zone, whose energy is on the order of the inverse lattice 
spacing . 

The evaluation of the loop integral over the gluon 4-momentum requires care in the 
choice of contour in the complex Euclidean energy plane. The path of the contour with respect 
to the poles in the propagators is determined by the coordinate space boundary conditions. The 
contour must be chosen so that the heavy quark propagator always vanishes for negative 
Euclidean time. Thus the propagator as a function of momentum and Euclidean time is given 



by 



dE 



2% v. 




iEa 




E v i sin a 



a i 



(18) 



1 



27U 



z v Q [z - 1] + T,v i s'mp i a 



(z = e iEa ) 
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X 




Figure 3 — Contour in the z = e plane for the one loop proper self mass 

and the necessity for this to vanish for negative t no matter the values of v and p implies that the 
contour in the z plane always encloses the pole. 

Although 3-momentum integrals will eventually have to be done numerically, it is 
necessary to do the energy integration analytically. With the notation 
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z = e 



iak A 



JL ~ k. a 
A = 4 Y, sin — 
i = l 2 



(lay 



B = Y v i s ' m (Pi + k i) a 
i = \ 



(19) 



c = E 



2 2 

V; COS 



i = l 



+ — 

v 2 y 



the Euclidean energy integral is 



2 2/p 4 a 
v Q e z + C 



271 ' z [- (z- 2 + l/z) + A][v (e' >4a z - 1) + B] 



(20) 



where the integration is around a closed contour in the z plane. Note that the quantities A and 
C are both positive, but B can take either sign. The singularities of the integrand are located at 



z + = 1 + — + - Ua + A 2 
± 2 2 



ip 4 a 



(21) 



1 - —B 



One of the pair of poles coming from the gluon propagator is within the unit circle in the z plane, 
and the other is outside. Depending on the sign of B, the pole coming from the heavy quark 



10 



propagator can be either inside or outside. The contour must be chosen to pass between the two 
gluon propagator poles, which is the ordinary proceedure. However, the requirement that the 
heavy quark propagator must vanish for all negative Euclidean times means that the the contour 
must always enclose the quark propagator pole, whatever the sign of B. 10 The appropriate 
contour is shown in Figure 3, for the non-standard case, negative B. 

The integrand decays sufficiently rapidly at infinity so that the energy integral is given 
by the residue of the single pole outside the contour, at z + . The resulting expression for the one 
loop proper self mass is 



E< 2 > =g 2 C 2 ^l 



(271 ) 3 



r 

2 2ip 4 a 

d 3 k 1? 1 (22) 

4A + A 2 [v (e lP4 \ - 1) + B] 



To facilitate comparison with the evaluation of the renormalization of the classical 
velocity by means of simulation, we express the perturbative result in terms of v- = v-/v , 

which can be seen from the reduced heavy quark Dirac equation (Eq. (15)) to be the actual 
transverse hopping expansion parameter of the lattice HQET. Its lowest order perturbative 
renormalization is given by 
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5„(ren) _ 1 
V- s V| . - V- = - _ 



(X. - v-.Aq) 



(23) 



Although this infrared finite expression can be numerically evaluated for any value of the 
classical velocity, explicitly evaluating the leading terms in a series on the components of v- 



brings out the structure of the result: 



3 2 

5v,. = c l v i + c 3 v. + c 12 v- Ev. + 



(24) 



The expansion coefficients are given by 



-s 2 c 2 

(27l) 3 



^ z + (cos9 ; . + z + - 2) 



4A + A 2 ( Z+ - l) 2 



(27l) 3 



3 (cos0. + 1)( Z+ - 3 cos9- + 2)(z + 2 + 2z + (cos0. - 2) + 1) 



2 \/4A + A 2 ( Z+ - l) 4 
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-g 2 c 2 

(27l) 3 



3 (cos0 J . + l)[2z + (z + -3cos0.+2)cos0 i . - (z + 2 Sz + + l)cos0. + z + ( z+ 2 -4 Z+ -3)] 



J J 



2 {AA + A 2 ( z+ - l) 4 



(j * i) 



(25) 



where the integration is over the periodic box - e [- n ,71 ]. 
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The numerical values of these coefficients, for g =6/(3 =6/6.1 are 

c l = - .23345568 

c 3 = - .04143250 (26) 
c n = - .03881061 

The difference between the values of the two cubic terms is a reflection of the reduced spatial 
symmetry of the lattice. The perturbative evaluation of the expansion coefficients of the 
corrections to the classical velocity are shown in Table II, for both naive and Lepage- 
Mackensie 11 tadpole improved mean field shifted values of the coupling. 

IV SIMULATION OF THE CLASSICAL VELOCITY RENORMALIZATION 

To evaluate the renormalization of the classical velocity non-perturbatively, we must be 
somewhat indirect. Without global gauge fixing, we cannot simulate the quark propagator. 
However, if we form a heavy-light meson, its classical velocity will be that same as that of its 
heavy quark component, and its propagator is easily computable in simulations. The calculation 
is greatly facilitated by noting that because the solution of the heavy quark lattice Dirac 
equation (15) is obtained by forward recursion, after n time steps the propagator is an (n - l) st 
order polynomial in the bare classical velocity components v- = v ; /v . 
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m, m m 



S iv) (t,p) = Y, v ~l ^2 v ~3 5(/,p,m) (27) 

OTj ,m 2 ,m 3 

The computation of the coefficients in this polynomial is highly efficient. The significance of 
the indices m ■ is that they are the maximum lattice displacement of the heavy quark propagator 

in the i th direction. The only essential difference between the asymptotic behavior of the heavy 
quark propagator and that of a heavy-light meson is the possibility, on the lattice, that the overall 
asymptotic normalization of the meson field composite operator depends on the residual 
momentum. 

s (v W) ~ c< v >(/0*-* (v)ao ' (28) 

The slope of the logarithmic derivative of the asymptotic propagator versus t, 
dS^\t,p)/d Pi \ p = dC^\p)ld Pi \ p = Q 3E^(p) 



S (v \t,p=0) C (v) (/? = 0) d Pi 

dC iv \p)/d Pi \ p = Q „ (phys) 



t 

p=0 



(29) 



C (v \p=0) 



v i 1 



is the negative of the physical classical velocity. Since actual lattices are finite, we use the lattice 
approximation to the momentum derivative 
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dS (v \t,p) 



AS (v \t,p) 



p=0 



p=0 



2Pr 



s (v) (t,p rPmin ) -s^>(t, Pr - Pmin ) 



(v), 



(30) 



Pmin 



2n 
~Wa 



The shift in the classical velocity to a given order n in the components of the bare classical 
velocity v t = v-/v , only depends on the heavy quark propagator at values of the spatial 

coordinates up to n lattice sites away from the initial heavy quark location. 



AS (v \t,p) 



A Pi 



p=0,O(v") 



Na 

471 \x,\<n 



2i sin . 



2 71 X; 



V 



N 



x S (v \t,x) 



(31) 



As with the perturbative evaluation, we display the results of the simulation to 3 rd order in the 



bare classical velocity. The coefficients of powers of v- ' in the logarithmic derivative are 
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AS (v '(t,p)/A Pi \ p __ 
S <v) (t,p=0) 



= E *, 



E *, 3 



AS(t,p,m i = 1)/Ap. |^ = 
S(f,p =0,m =0) 

AS(t,p,m i = 3)/Ap i |- = 
5(f,p =0,m =0) 

AS(t,p, mi = l)/A Pi \ p = Q S(t,p=0, mi =2) 



S(t,p=0,m=0) 2 

AS(t,p,m i = 1 ,m ; - =2)/Ap- |^ = 
S(t,p=0,m=0) 

AS(t, P , mi = 1)/Ap. - = S(t,? = 0,mj=2) 
S(t,p=0,m=0) 2 



(32) 



Averaging over directions of course improves the signal to noise ratio. 

We have simulated these propagator expansion coefficients on an ensemble of lattices and 
Wilson light quark propagators made available to us by the Fermilab ACP-MAPS 
Collaboration . These consisted of 98 lattices of size 24 x 48 with lattice coupling (3 = 6.1 
along with Wilson quark propagators with hopping constant K = .154. We computed the 
coefficients of the leading terms in the v ■ expansion of the heavy quark propagator, and formed 

heavy-light meson propagators (or, more precisely, their v- expansion coefficients) in the 

standard manner. We did not use smearing or otherwise improve the saturation of the resulting 
objects by the lowest mass state contribution, but it is our intention to use systematic variational 
methods to improve the isolation of the ground state in a future, more refined analysis ' . 
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The simulation did take advantage of the fact that, since the direction of time is 
conventional, purely forward, or purely backward propagation of the heavy quarks in Euclidean 
time gave effectively independent ensembles of lattices and propagators. 

Figure 4 shows the results of the simulation of the heavy-light meson propagator. The 
solid circles are the values of the leading, zeroth order term in the expansion of the propagator 
in powers of v t , which is the static propagator. Also shown are a single exponential fit to the 

static propagator with a fall-off given by the average efffective mass over the range between 2 
to 6 units of Euclidean time, and the fractional saturation of the propagator by this single 
exponential. The other propagator components shown are those coefficients of the linear and 

cubic terms in the the v- expansion of the propagator that enter into the simulation of v^ p ys ^ 

through third order in the bare classical velocity components, that are shown in Eq. (32). Table I 
shows the fitted slope of each coefficient. 

The degree to which a single exponential saturates the propagator is an indication of how 
well the ground state is isolated. The fact that there good saturation for at most a very small 
range of Euclidean times shows that by the time one has reached sufficiently large Euclidean 
separations to isolate the ground state, the statistical errors in the propagator have become 
overwhelming. This is an independent way of seeing that there is not an extended region of 
Euclidean time over which the effective mass of heavy-light propagator is constant 14 . A plot of 
the effective mass of the m = propagator in each Euclidean time interval is shown in 
Figure 5. The absence of a stable plateau must be remediated either by better statistics, or 
probably more effectively by improving the heavy-light meson wave function. 
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Figure 4 — The coefficients of powers of the bare classical velocity in the heavy-light 

meson propagator through third order 



The three combinations of ratios of propagator components that enter into the third order 
evaluation of the physical classical velocity shown in Eq. (32) are shown in Figure 6. The errors 
are the propagated statistical errors associated with the individual terms in the propagator 

18 



1.3 
1.2 
00 1.1 



CO 
CO 



> 10 

o 
CD 

3= 0.9 

LU 



0.8 
0.7 



_ • 

• 

• 

■ 




i 

f i - 

< 

i i i i i i i 


► < 

- i 


> 

i - 


- i i 



23456789 
Euclidean Time 



10 11 12 



Figure 5 — The effective mass of the zero classical velocity, zero residual momentum 
heavy-light meson propagator between each time slice 



expansion. As shown by Eqs. (29) and (32), the asymptotic slopes of these combinations of 



propagator contributions are the coefficients of the indicated terms in the expansion of - v ■ 



(phys) 



in powers of v . 1 . In order to estimate them, we fit each coefficient function to a straight line, 

starting at t = 3 for the linear coefficient and at t = 4 for the cubic coefficients. The fitted slopes 
are shown in Table I. 



V TADPOLE IMPROVEMENT 

In the analysis so far, we have not taken advantage of the tadpole improvement 
proceedure of Lepage and Mackensie 11 . This consists in simply replacing each link variable 
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Figure 6 — The three combinations of propagators that are the coefficients of the linear and 
cubic terms (in v) in the logarithmic derrivative of S^ v \ as in Eq. (32) 



U (n) by that link divided by its average value U (n) I u Q , where the average value of the link 



is given by the gauge-invariant expression 
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Table I - Fitted slopes from simulation of heavy-light meson propagator 



Coefficient of 


Fitted Slope 


Fitting Range 


v i 


-.713 + .092 


3 - 9 




.982 + .087 


4 - 8 


v f v/ (i * j) 


.764 + .078 


4 - 8 



u 



-{TrU) 
3 



1/4 

(33) 



where □ denotes the product of the link variables around an elementary lx 1 plaquette. In the 
present context there are two aspects to this proceedure. One is the identification of the 
appropriate perturbative coupling associated with a simulation at a given value of (3, which would 
simply be 6/g without tadpole improvement. For the Wilson action used to generate the 
Fermilab lattices, 



4^4 rrD ^"^r^T rrn (34) 

8 □ 3 8 Uq □ 3 



rs A 

(3 instead is identified with 6/ g u . Since u Q is necessarily less than one, this has the effect 

of increasing the value of the effective perturbative coupling corresponding to a lattice simulation 
at a given value of p\ 

The other aspect is the improvement of the heavy quark propagator. If each link in the 
discretized reduced Dirac equation (15) is replaced by itself divided by w , then modulo an 
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overall normalization, which of course does not affect the determination of the physical classical 
velocity, the effect of all the factors of u Q can be absorbed into a redefinition of the lattice bare 

classical velocity and a mass renormalization. Specifically, if we absorb a factor of 

— yi 

Uq ' = exp (- n t \uuq) into the heavy quark propagator, we recover the original discretized 
equation (Eq. (15)) for the modified propagator function, but with the replacement 

v • v ■ 

v t = _L => _L (35) 

-fl 

The factor of u Q ' is evidently just a shift in the mass by - lnw . As with any other residual 

mass, it is without physical effect. The replacement in the classical velocity means that it is 
actually the quantity v • / u Q which is the transverse hopping constant. The factors of w that 

come from replacing v • by v • / u Q can be absorbed into the coefficients of powers of the bare 

classical velocity v- in the expansion of the physical classical velocity v^ p ys \ The coefficient 

of a term with N powers of v- , whether determined by simulation or perturbatively, is 

multiplied by u Q . In terms of the expansion coefficients of 8 v of Eq. (24), if we denote the 

improved coefficients by c, the relation between the improved and unimproved sets of 
coefficients is 
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Ci = 



+ 1 - Uq 



C-i = 



c 12 " 



12 



(36) 



i+j + k 



(i + j + k > 1) 



Table II — Comparison of Simulated and Perturbative Expansion Coefficients 



Coeffi- 
cient 


Term 


Naive Coupling 


Tadpole Improved 


Simulation 


One-Loop 


Simulation 


One-Loop 


c l 


v i 


-.287 + .092 


-.2334 


-.186 + .105 


-.3113 


c 3 




-.982 + .087 


-.0414 


-1.46 + .130 


-.1048 


c 12 


v,v/ 

0" * i) 


-.764 + .078 


-.0388 


-1.14 + .116 


-.0982 



V DISCUSSION OF THE RESULTS 

It is striking to compare the perturbative evaluation of these expansion coefficients with 
the results of the simulation performed at (3 = 6.1 . This is shown in Table II, both with and 
without tadpole improvement. The coefficients of the linear term agree rather well between the 
simulated and perturbative evaluations. There is a statistically insignificant preference for the 
use of the naive rather than the tadpole improved coupling in this context. However, the 
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Figure 7 — Composition of the v- term in the heavy-light meson propagator 



perturbative and simulated cubic coefficients are in complete disagreement. It is difficult to 
surmise which estimate, if either, is correct. While it is possible that the lattice coupling is 
simply too large to neglect higher order terms, is should be pointed out that the simulations of 
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Figure 8 — Composition of the v- v- term in the heavy-light meson propagator 



the cubic (and higher) coefficients are very delicate calculations. 

Each of the cubic coefficients in the expansion of the logarithmic derivative of the heavy- 
light meson propagator is the difference of two terms, each of which is substantially larger than 
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their difference. Furthermore, each of the two terms has a P asymptotic behavior, while their 
difference should only grow linearly. Thus as one tries to eliminate non-asymptotic artifacts by 
going to large t, the statistical quality of the calculation rapidly breaks down. Fortunately, the 
dominant, linear coefficient is free of this particular difficulty. The large cancellations are 
roughly of the same magnitude in both cubic coefficients. The degree of cancelation is shown 
for the v- 3 term in Figure 7 and for the v-Vy 2 term in Figure 8. Note that the P asymtotic 
behavior of the individual terms is clearly visible from the smallest values of the Euclidean time. 

Let us sum up the results of these analyses of the renormalization of the classical velocity 
in lattice version of the heavy quark effective theory. The origin of this new renormalization is 
the reduction of the Lorentz invariance of the continuum to hypercubic symmetry on the lattice. 
Both in perturbation theory and via simulations the shift in the classical velocity is quite 
substantial. The most reliable results concern the first-order, multiplicative shift. The two 
analyses are in fairly good agreement on this point. Since the lattice coupling corresponds to a 
continuum coupling constant near 1, it is probably most reasonable to use the result coming from 

the simulation, which is intrinsically non-perturbative. This first-order (in v ■ ) multiplicative 
renormalization is 



„(phys) „(bare) 

„ (bare ) 
V i 



— < 



.186 + .105 {Tadpole improved) 

.287 + .092 (Unimproved) 



(37) 
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In contrast to the first-order situation, the simulated and the perturbative third-order results 
differ from each other by an order of magnitude. Neither seems to us terribly reliable, though 
it is hard to invent a reason for the perturbative result to be an order of magnitude off. By 
contrast, the substantial intrinsic errors we have described above in the simulation of the cubic 
coefficients make it quite conceivable that the asymptotic slope has not yet been observed in the 
simulations reported here. Since the perturbative result for these coefficients is tiny, one may 
hope that for moderate values of the classical velocity, the first-order correction suffices for 
practical use. 
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